////////////////////////////TCF/////////////////////////////
// TPC Coherent Fit (TCF)                                 //
// Author: Xianguo Lu, lu.xianguo@gmail.com               //
// Reference: Xianguo Lu, "Exploring the performance      //
// limits of the ALICE Time Projection Chamber and        //
// Transition Radiation Detector for measuring identified //
// hadron production at the LHC", dissertation at the     //
// University of Heidelberg (2013).                       //
////////////////////////////////////////////////////////////

//#include "headers.hh"
#include "style.h"

//drawing outmn2d_post_draw2d.root
//in case of outmn2d_pre_draw2d.root, change every post to pre

/*
//good = to be drawn

  KEY: TH2D     outmn2d_post_original;1 horiginal tpcsig : log10p (inequal x bin width) //good
  KEY: TH2D     outmn2d_post_fit;1      horiginal tpcsig : log10p (inequal x bin width)
  KEY: TH1D     outmn2d_post_mean0;1 //OK, but not nice to draw
  KEY: TH1D     outmn2d_post_mean1;1 //OK, but not nice to draw
  KEY: TH1D     outmn2d_post_mean2;1 //OK, but not nice to draw
  KEY: TH1D     outmn2d_post_mean3;1 //OK, but not nice to draw
  KEY: TH1D     outmn2d_post_chi2;1 //OK, but not nice to draw
  KEY: TH1D     hfitdatax170;1  x -1.525000e-01 //good
  KEY: TH1D     hfitdatax170fit;1 //good
  KEY: TH1D     hfitdatax170f0;1 //good
  KEY: TH1D     hfitdatax170f1;1 //good
  KEY: TH1D     hfitdatax170f2;1 //good
  KEY: TH1D     hfitdatax170f3;1 //good
...
  KEY: TH1D     hfitdatax334;1  x 1.275000e+00 //good
  KEY: TH1D     hfitdatax334fit;1 //good
  KEY: TH1D     hfitdatax334f0;1 //good
  KEY: TH1D     hfitdatax334f1;1 //good
  KEY: TH1D     hfitdatax334f2;1 //good
  KEY: TH1D     hfitdatax334f3;1 //good
  KEY: TGraph   outmn2d_post_mean0gr;1 //good
  KEY: TGraph   outmn2d_post_mean1gr;1 //good
  KEY: TGraph   outmn2d_post_mean2gr;1 //good
  KEY: TGraph   outmn2d_post_mean3gr;1 //good
  KEY: TGraph   outmn2d_post_chi2gr;1 //good
 */

TString fTag;
TString fDir;

void singledraw(const Int_t ix, const Int_t ntype)
{
  //============================================================================================
  //                          measured distribution
  //============================================================================================
  TH1D* hraw=(TH1D*)gDirectory->Get(Form("hfitdatax%d",ix));
  if(!hraw){
    return;
  }

  TCanvas *c1=new TCanvas(Form("c1%d",ix),"",600,400);
  gPad->SetLogx(0);
  
  //xx: log10(p)
  TString tit(hraw->GetTitle());
  tit=tit(tit.First(' ')+1, tit.Length());
  if(tit.Length()==0){
    printf("bad tit!!\n");    exit(1);
  }
  const Double_t xx=tit.Atof();
  printf("hraw: %s, tit: %s, xx: %f\n", hraw->GetTitle(), tit.Data(), xx);

  gPad->SetLogy(xx<0);
  if(gPad->GetLogy()){
    hraw->SetMaximum(1);
  }

  //range in dE/dx, not of log10 p!!
  hraw->GetXaxis()->SetRangeUser(30, 100);
  if( TMath::Power(10,xx) < 1.0 ){
    hraw->GetXaxis()->SetRangeUser(20,160);
  }
  hraw->SetMarkerStyle(24);
  hraw->SetMarkerSize(1);
  hraw->SetLineColor(kBlack);
  hraw->SetXTitle("TPC signal (arb. units)");
  hraw->SetYTitle("Entries");

  hraw->Draw();

  //============================================================================================
  //                          individual particle fit
  //============================================================================================

  TLegend * lg=new TLegend(0.66, 0.58, 0.94, 0.88);
  lg->SetFillStyle(0);
  lg->SetHeader(Form("p=%.2f GeV/c", TMath::Power(10,xx)));

  const Int_t col[]={kBlue, kRed, kGreen+3, kMagenta, kYellow};
  const Int_t fillsty[]={3004, 3005, 3006, 3007};
  const TString types[]={"p",  "#pi", "K", "e"};
  for(Int_t itype=0; itype<ntype; itype++){
    TH1D* hfitii=(TH1D*)gDirectory->Get(Form("hfitdatax%df%d",ix, itype));
    if(!hfitii){
      printf("no hfitii!! %d %d\n", ix, itype); exit(1);
    }
    hfitii->SetFillColor(col[itype]);
    hfitii->SetFillStyle(fillsty[itype]);
    hfitii->SetLineColor(col[itype]);
    hfitii->Draw("same");
    lg->AddEntry(hfitii, types[itype],"lf");
  }

  //============================================================================================
  //                          total fit
  //============================================================================================

  TH1D*hfitall=(TH1D*)gDirectory->Get(Form("hfitdatax%dfit",ix));
  if(!hfitall){
    printf("no hfitall!! %d\n", ix); exit(1);
  }
  hfitall->SetLineWidth(2);
  hfitall->SetLineColor(kBlack);

  hfitall->Draw("same");lg->AddEntry(hfitall,"fit","l");
  lg->Draw("same");

  c1->Print(Form("%s/%s_slice%d.eps", fDir.Data(),  fTag.Data(), ix));
}

void plotdraw2d(const TString fin="outmn2d_post_draw2d.root", const Int_t ntype=4)
{
  const TString tmp=fin(0,fin.First('.'));
  fDir=Form("eps_%s", tmp.Data());
  gSystem->Exec(Form("mkdir -p %s/", fDir.Data()));

  fTag = fin.Contains("pre")?"pre":"post";

  TFile::Open(fin);
  //========================================================================================================
  //                       draw 2D measured TPC vs p distribution overlaid with fitted mean from particles
  //                       outmn2d_post_original + outmn2d_post_mean%dgr 
  //========================================================================================================
  TCanvas* c2=new TCanvas("c2","",600,400);
  gPad->SetLogx();
  gPad->SetLogz();

  TH2D* hfit=(TH2D*)gDirectory->Get("horiginal");
  ToNaturalScale(hfit);

  hfit->SetXTitle("p (GeV/c)");
  hfit->SetYTitle("TPC signal (arb. units)");
  hfit->SetTitle("Entries");
  hfit->GetXaxis()->SetRangeUser(0.2, 40);
  //hfit->SetMinimum(1e-5);
  //hfit->SetMaximum(1e-1);
  
  hfit->Draw("colz");
  
  const Int_t linesty[]={kSolid, kDotted, kDashed, kDashDotted};
  //const Int_t col[]={kBlue, kRed, kGreen+3, kMagenta, kYellow};
  
  for(Int_t itype=0; itype<ntype; itype++){
    TGraph* hmean = (TGraph*)gDirectory->Get(Form("outmn2d_%s_mean%dgr", fTag.Data(), itype));
    ToNaturalScale(hmean);

    hmean->SetLineWidth(2);
    hmean->SetLineColor(kBlack);//col[itype]);
    hmean->SetLineStyle(linesty[itype]);
    hmean->Draw("pl");
  }

  c2->Print(Form("%s/%s_dist2d.eps",fDir.Data(), fTag.Data()));

  //========================================================================================================
  //                       draw reduced chi2 vs p
  //                       outmn2d_post_chi2gr
  //========================================================================================================

  TCanvas* c4=new TCanvas("c4","",600,400);
  gPad->SetLogx();

  TGraph* chidof = (TGraph*)gDirectory->Get(Form("outmn2d_%s_chi2gr", fTag.Data()));
  ToNaturalScale(chidof);

  chidof->SetMarkerStyle(24);
  chidof->SetMarkerSize(1);
  chidof->SetMarkerColor(kBlack);
  chidof->SetLineColor(kBlack);
  chidof->GetXaxis()->SetTitle("p (GeV/c)");
  chidof->GetYaxis()->SetTitle("#chi^{2}/NDOF");
  
  chidof->SetMinimum(0.0);
  chidof->SetMaximum(3);
  chidof->Draw("apl");

  TF1 * f1=new TF1("f1","1",-1000,300);
  f1->SetLineColor(kBlack);
  f1->SetLineWidth(1);
  f1->Draw("same");
  
  c4->Print(Form("%s/%s_chi.eps", fDir.Data(),  fTag.Data()));

  //========================================================================================================
  //                       draw single slices TPC distribution with fit signal shapes
  //                       hfitdatax: measured, total fit, individual particles
  //========================================================================================================
  for(Int_t ii=1; ii<=hfit->GetNbinsX(); ii++){
    singledraw(ii, ntype);
  }
}

int main()
{
  return 0;
}
